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DISCRETE CALCULUS OF VARIATIONS FOR QUADRATIC 
LAGRANGIANS. CONVERGENCE ISSUES 

PHILIPPE RYCKELYNCK* LAURENT SMOCH* 

Abstract. We study in this paper the continuous and discrete Euler-Lagrangc equations aris- 
ing from a quadratic lagrangian. Those equations may be thought as numerical schemes and may 
be solved through a matrix based framework. When the lagrangian is time-independent, we can 
solve both continuous and discrete Euler-Lagrange equations under convenient oscillatory and non- 
resonance properties. The convergence of the solutions is also investigated. In the simplest case of 
the harmonic oscillator, unconditional convergence does not hold, we give results and experiments 
in this direction. 
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1. Introduction. The principle of least action may be extended to the case 
of non-diffcrcntiable dynamical variables by replacing in the lagrangian £(x, x) the 
derivative x(i) of the dynamical variable x(t) with a 2N + 1-terms scale derivative 

N 

□«x(t) = (Hx(t + ie)x-i(t), te[a,b], (1.1) 

i=-N 

see O IH |6] . Here, e stands for some time delay and Xi CO denotes the characteristic 
function of the interval [max(a, a+ie) , min(&, b+ie)] . Critical points of classical actions 
are characterized by the classical Euler-Lagrange equations V x £ — d/diV x £ = 0. 
Similarly, we proved in [6] that the equations of motion for discretized actions are 

V x £ + D-eVx/: = 0. (1.2) 

We abbreviate as C.E.L. and D.E.L. the classical and discrete Euler-Lagrange systems 
of equations respectively. 

In this paper we work with lagrangians of the shape £(x, x) and £(x, D E x) where 
£ : C d x C d — > C is a quadratic polynomial. We are interested in solving C.E.L. and 
D.E.L. under Dirichlet conditions. More accurately, we study the existence and the 
unicity of pseudo-periodic solutions z(t) of C.E.L. and y e (t) of D.E.L., e being fixed. 
The underlying assumptions for this to occur may be thought as an "oscillatory" con- 
dition for the lagrangian £(x, y) and as a "non-resonance" condition for the Dirichlet 
problem associated to C.E.L. and D.E.L.. With this in mind, we address the problem 
of convergence of y E (i) to z(t). 

The paper is organized as follows. Section 2 gives notation and basic definitions 
used throughout. In Section 3, we develop a matricial based framework to solve 
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D.E.L. for all quadratic time-dependent lagrangian. In Section 4, we provide under 
mild assumptions formulas for the components of the pseudo-periodic solutions of 
C.E.L. and D.E.L. when the lagrangian does not explicitly depend on time. This 
allows us to compute in some particular cases the phases of y e (i) and z(t) help to 
the matrix A. Section 5 is a preliminary discussion of convergence of y e (t), uniformly 
locally in ]a,b[, as e tends for stationary lagrangians and pseudo-periodic solutions. 
If £ is a non-resonant oscillatory lagrangian and D e is a well-chosen three-terms 
operator, the previous convergence property is the content of our main theorem which 
is proved in Section 6. In Section 7, we give numerical experiments to illustrate the 
non-unconditional convergence of solutions. 

2. Preliminaries. First, let us collect some notation and definitions from [6]. 
If the context is clear enough, i = \f—l. Let [a, b] be some interval of time and a 
time delay e > be fixed throughout. The integers d and N denote respectively the 
"physical" dimension and the number of samples in C d . We define for to € [a, 6] the 
grid Gt ,e — {to + n£N}n [a, b]. We denote by Id the identity matrix of size d. 

Let C pw (d,N) be the space of the functions x : [a, b] C d continuous on each 
interval [a + ie, a + (i + l)e] n [a, b] for all i £ {~N, . . . , N}. The two functional spaces 
C 1 ([a, b], C d ) and C pw (d,N) are Banach algebras with uniform norms. The operator 
□ e given in (jl.ip is a continuous linear endomorphism of C pw . 

Now, let be given six mappings P,Q,J X € C x {[a, 6], C dxd ), J 2 , J 3 G C l ([a,b],C d ) 
and J4 £ C 1 ([a, b],C). We suppose that for all t, P(t) and Q{t) are symmetric and 
J \{t) is skew-symmetric. We set 



and we define the quadratic lagrangians £(x,x) and £(x,D e x). If the coefficients in 
(|2.1j) do not depend explicitly on time, we shall say that C is stationary. 
We will consider actions *4 C ont(x) and Adisci'x) of the shape 



The actions A cont ■ C 1 ([a, 6], C d ) —> C and Adisc ■ C pw (d, N) —> C are continuous and 
Frechet differcntiable everywhere. 

We give in [6] the necessary first order conditions of local optimum of A con t and 
Adisc under the Dirichlet constraints x(a) = d a and x(6) = df, in the previous spaces, 
where d a and db are two fixed vectors in C d . The Euler-Lagrange equations associated 
to each action in (|2.2p can be written as 



£(x, y) = - Vy + - T xQx + T xJ iy + T J 2 y + T J 3 x + J 4 . 



(2.1) 





(2.2) 



- Px + (-P + 2Ji)x + (Ji + Q)x - J 2 + J3 = 0, 
□_ £ (PD £ x) - □_ s (Jix) + JiD e x + Qx + D_ e J 2 + J 3 = 0. 
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(2.3) 
(2.4) 



The problem of convergence as e tends to of the operator in the l.h.s. of (|2.4j) 
to the corresponding operator in the l.h.s. of (|2.3I) has been studied in [6]. In this 
context, we introduced the class of discretization operators given by 

□ M X (i) = -X+i(*)|x(t - e) + ^x(t) + X -i{t)-x{t + e). (2.5) 

where r,s £ C. In fact, (|2.5p gives the shape of three-terms operators satisfying 
□ e (l)(i) = inside [a + 2e, 6 - 2e]. 

Without assuming the convergence of the schemes in the previous sense, we focus 
on the following two problems. Are the Dirichlet problems for (j2.3[) and (|2.4I) well- 
posed? Are there periodic or pseudo-periodic solutions? In fact, if e = (b — a)/M for 
some M £ N* and if x e (t) is a solution of D.E.L., then x e is uniquely determined on 
the grid Q a ,e- We shall see later how to construct from x e the unique corresponding 
pseudo-periodic solution y e of (|2.4I) . 

3. An effective method for solving D.E.L.. In this section, the datas N, d, C, 
d a , df,, £, D e are fixed but arbitrary. 

3.1. D.E.L. as delayed functional equations. The equations (|2.4|) may be 

thought as a mixture between recurrence equations and delayed functional equations. 
Let us transform the problem of solving (|2.4j) into an infinite set of problems, each 
of one dealing with recurrence vector equations with the additional difficulty of the 
perturbation of the boundaries. We can identify each function x : [a, b] — > C d to the 
infinite set of finite sequences (x(io + ne), n £ Z) with indices such that ^^r 2 - < n < 
and where to lies in an interval of length e. In fact, because (|2.4[) involves the 
second order operator □_ e D e , it may be formulated in an abstract manner as 

F(i, x(f - 2Ne) 1 x{t - (2N - l)e), . . . ,x(i + (2N - l)s), x(f + 2Ne)) = (3.1) 

where F : C 1+(4Ar+1)d -> C d contains the coordinates of the l.h.s. of (|23]l . Hence, 
we solve (|3.ip with respect to x(to + 2Ne) for fixed io, or what amounts to the same 
thing by expressing x(t + 4Ne) as a function of x(i + he) for k £ {0, . . . , 4iV — 1}. 

For instance, the case N = 1 and d arbitrary is the most interesting one, and we 
may rewrite in this case the equations (|2.4p as a system of d equations 

Xi(t-2e) [ciC-iXiX2Pij(t-e)] 
+Xi(t - e) [coCiXiPij(f - e) + coc-ixoXi-Rj W + cixi^ijjif* - e) + c-iXi(Ji)«(*)] 
+!Ci(t) [c?XoXi-P S j(i - e) + c§xoJ\i(t) + c 2 _ lXoX _ 1 P, J (f + e) + Xo<5ij(*)] 

+ e) jcocixoX-i-Pw(t) + c c_iX-iPij(f + e) + ciX-i(Ji)«(t) + c-iX-i(Ji)ji(t + e)] 
+a;i(t + 2e) [c-iciX-iX-2-Py (t + e)] 

+(ciXoXi(Ja)j(* - e) + coXo^jM + c_iXoX-i(J2)j(* + e) + Xo(Ja)i(*)) = 

(3.2) 

for each j £ {1, . . . , d}, Vt G [a, 6], with summation on i when repeated. This equation 
has been heavily used for numerical experiments. 
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3.2. Solving D.E.L. in the safety interval. Given i n € [a, 6], we define the 
safety interval as the segment Is C N such that 



n e l s iff h + (n - j)e € [a, b] for all j G {0, . . . , 4JV — 1}. (3.3) 
We convert now (|2.4|) into a linear recurrence in <C 4dN . For n € Z,g , we set 

/ x(t + ns) \ 
x(t + (n - 

e C 4dJV . 

Vx(t + (n-4iV+l)e)/ 

When n E Is, every characteristic function occuring in (|2.4j) equals to 1. Then, there 
exists well-defined matrices A n € (Q4d.ivx4d.zv anc j vec t ors Dn g C 4dAr , depending only 
on n, D e , P, Q, J\, J2 and J3, such that (|2.4[) is equivalent to 

v ?!+ i = i„v„ + b„. (3.4) 

The matrix A n is defined at this stage if n,n + 1 € I s , and admits a block structure 
with 4iV x 4iV blocks of size d x d. On block rows 2, 3, . . . , N, the blocks are either 
identity blocks or zero blocks, and on block row 1, the blocks Bi n: i G {1, . . . ,4iV} 
will express the matricial coefficients in the equation derived from (|2.4p by solving it 
w.r.t. x(ta + ne). In this way, the matrix A n is the block companion matrix of the 
matrix polynomial 



T YiN o v4N-1 d y4JV-2 r v TJ 

— i>l,nA — -t>2,n^ — ... — ±>4jv_i !rl A — r>AN,n- 

For sake of clarity, if C is stationary, it turns out that those AN blocks have the 
shape 



Bi, n = c'ih + c'lP- 1 ^ + c'!'P- 1 Q 



(3.5) 



where the constants c^,c",c"' depend only on i and the coefficients Cj. Moreover, if 
N = 1, the following formulas for T A n and ^hn display the general structures of A n 
and b„ 



T A r . 



I, 



c i+c_i)c n 
cic_i 

cic_i 
(ci+c-i)co t 



/' •/, h 

P- X Q I d 



-/ d 

(-£r(^n- £ (i) + J 3 ) 0) 



cic_i 

L_ D i / 

cic_i 





Id 
0/ 



(3.6) 



(3.7) 



3.3. Conditions for D.E.L. to be well-posed. Let us consider the problem 
of solving D.E.L. under Dirichlet conditions. In the following result, we deal with 
existence, uniqueness and determination of the restrictions of the solutions of D.E.L. 
to the various grids Gt ,e- 

Theorem 3.1. Let t e [a, b] and e > 0. 
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• If {a, b} C Gt ,e> either there does not exist any solution x e (t) on [a, b], or the 
restriction of each solution to Gt ,e is uniquely determined by the vectors d a 
and d;, in C d . 

• If for instance {a, b} D Gt ,s = {&}, then the set of solutions x £ : Gt ,s >-> C d 
of \2.J$ is in one-to-one correspondance with C d . 

• If {a, b} PiGto.s = 0, £/ien i/ie set of solutions of \2.1$ on Gt ,e is in one-to-one 
correspondance with C d x C d . 

Proof. We assume that N = 1 only to be more explicit, the case iV > 1 having 
the same qualitative features. Let us suppose that {a, b} C Gt ,s and w.l.o.g. that 
to = a and b — a = Me where M £ N*. For the need of the proof, we pursue the 
construction of A n when n $.%s- In that case, some characteristic functions occuring 
in (|2.4p vanish, this relationship is no more of order d, and the sizes of v ra and A n 
must change. We have x e (a) = d a and we set x e (a + e) = d s e C d which is introduced 
without being determined at this stage, firmly from recurrences. Plugging t — a in 
recurrence p. 21) and solving, we first get 

x E (a + 2e) = B hl d s + B 2A d a 
where Bii,B 2 i are blocks similar to those occuring in (|3.5[) . Next, with £ = a + e we 
find 

x E (a + 3e) = (B h2 B hl + B 2a )d s + (B h2 B 2A + B 3 , 2 )d a . 
The following iterations express x e (a + ne) as a linear combination of the vectors 
d s , d a , with coefficients being polynomial matrices in Bi^. From index from n = 4 to 
n = M — 3 the recurrence (|3.2p becomes or order 4iV + 1 and may be reformulated as 
(|3.4p . Finally, the three last steps n = M — 2, M — 1, M are similar and imply three 
systems of decreasing sizes. In order to convert matricially this process, we introduce 
the five rectangular matrices Ai 



B i,i 


B% 


Id 








Id 



Ai=\ I d \,A 2 



( B iy2 B 2y2 B 3 , 2 \ 

Id 

o i d o 

V o o / d y 



Bl,M-3 B2.M-3 -B3.A/-3 — -f<2 
A M -3 = I Id I , A M -2 

Id 



Am-i — {B\ — Id) ■ 



B\ t M-2 B 2 _M-2 —Id 

Id 



The operators Ai, A 2 are used to compute the values x e (a + ne) for n — 2, 3 linearly 
as functions of d a , d s . Next, we have 

x £ (a + ne) = (A„-i ■ • ■ A 3 )A 2 A, ( J j (3.8) 

for 4 < n < M — 3. Finally, Am-3, A/w-2> Am-i are used to find x e (a + ne) for 
M — 2 < n < A/. At the end of the process, we get the shooting equation for the 
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vector d s : 



d 6 = x e (b) = x e (a + Me) = A M -iA M -2A M -3( a m-4 ■ ■ ■ A Z )A 2 A 1 ( ^ s J . (3.9) 

Now, existence and unicity of the restriction of x £ to the grid Q a £ = Gb,e is equivalent 
to the fact that the shooting method is successful, that is 

det (am^Am-2 ■ • • A 2 A 1 ( )) ^ °' (3 ' 10) 

Let us consider now the cases where (b — a)/e is not an integer so that \Gt .e H 
{a, b}\ < 2, the previous matrix formalism being similar. If {a, 6} D Gt ,e = then 
any vector d s £ C d determines a solution x £ (i) on Gt ,e- The case {a, b} n Gt ,e = {b} 
is entirely similar and we have infinitely many choices for d s = x e (6 — e). Lastly, if 
{a, 6} n Gt ,e = then we first may choose arbitrarily the two vectors x e (min Gt ,e) 
and x e (e + mmGt ,s) m C d and we use (|2.4p to compute iteratively the values of x E 
on Gt ,e- ■ 



Remark 3.1. Let us note that if £,a, 6, are fixed, the underlying determinant of 
Am-a ■ ■ ■ A% is a nonzero polynomial of degree less than 2d x (4M — 6) w.r.t. the 
coefficients c_i,co,ci and does not vanish generically. 

3.4. Eigenvectors of the matrix A n when n G Is- As it is the case for the 
sequences of vectors satisfying ordinary linear recurrences, the qualitative features of 
the solution x e (< + ne) of D.E.L. are reflected by properties of the spectrum Sp(A n ) 

Of A n . 

PROPOSITION 3.2. The eigenvectors of A n in C idN have the shape 
V- T (wA 4Ar -\wA 47V - 2 ,...,w) ; 

4N 

where w £ ker(J^ B^X^^ - X iN I d ) C C d . We have det(A n ) = {-l) d and 
i=i 

AN 

det(A n - \I idN ) = det(^] B hn X 4N - 1 - X iN I d ). 

i=l 

Proof. The two results are well known in the scalar case d = 1. Let us give some 
details when we deal with characteristic functions and d > 1. 

If v e c 4 dAf j g an eigenvector of A n associated to A € C, we partition it as T v = 
T (w4N, . . . , wi) where w; E C d . We next identify the corresponding blocks of size 
d x 1 in A n \r = Av to get w t = Awj_i = A !_1 wi for 2 < i < AN. Renaming wi as w 
and plugging the vectors in the first block row of A n v yield the first property. 
The second property may be easily proved by using matricial techniques for parti- 
tioned matrices (see for instance [71 pp. 36]). ■ 
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4. Pseudo-periodic solutions of C.E.L. and D.E.L. for stationary la- 
grangians. In this section, the datas N, d, d a , dj, , e are fixed but arbitrary, and C 
is stationary. We will say that £ is a stationary non-resonant oscillatory lagrangian 
w.r.t. the datas N, d, d a , db, e, O s if and only if D.E.L. and C.E.L. admit one and only 
one pseudo-periodic solution y e (t) and z{t) respectively. 

4.1. Solving C.E.L.. Let us study first the existence, unicity and periodicity 
or pseudo-periodicity of the solutions of ()2.3|) . 

Proposition 4.1. Suppose that C is stationary and that for some matrices 
n u Q 2 & C dxd we have 

Pn 2 + 2iJ x n + Q = and det(exp (i(b - a)fi 2 ) - exp (i(b - a)Qi)) ^ 0. (4.1) 

Then, for all d a ,db G C d there exists one and only one solution of C.E.L. \2. 3\) to- 
gether with Dirichlet boundary conditions. Moreover, iffl\ and f^2 are diagonalizable, 
each component f(t) of z(t) may be written as 

f(t) = cst+ ^ X! cstkj,u(dk)jexp(icj(t- a)). (4.2) 

k = a,b weSp(0!)uSp(0 2 ) 
3 6 [M 

where the various constants depend only on their indices as well as b — a and the 
eigenvalues of tti and f2 2 . 
Proof. We see first that 

z(t) = exp (i(t — a)f2i)zi + exp (i(t — a)J7 2 )z 2 — Q^ 1 J3, (4-3) 

is a solution of (|2.3j) for all zi,z 2 € C d . In order to fit the Dirichlet conditions, 
the vectors Zi,z 2 must satisfy Zi + z 2 = d a + Q~ 1 J^ and exp (i(b — a)Qi)zi + 
exp(i(b — a)fi 2 )z 2 = db + Q^ 1 ^. Due to f|4. 1 [) . the previous system is Cramer and 
the solution is equal to zi = i?e 2 and z 2 = —Rei where R € C dxd and ei,e 2 £ C d 
are respectively defined by 

R = (exp(i(6 — a)0 2 ) — cxp(i(b — a)^))^ 1 and 
e fc = exp(i(b - a)fl k )d a - d b + (exp(i(b - a)Q k ) - 7 d )Q _1 J 3 . 
By considering the previous formulas, we see that each component f(t) of z(t) depends 
linearly on (d a , db, Q _1 J3) g C 3d and may be returned as (|4.3[) where the constants do 
not depend on t, d a , db nor on Q^ 1 J3. Indeed, since fifc is diagonalizable for k = 1, 2, 
each entry in exp(iif2fc) is a monomial exponential w.r.t. t. Thus, each component of 
(j4~3| has the shape (|L2|) . ■ 

As (|4.3I) shows, the solution z(t) of C.E.L. is pseudo-periodic if and only if the 
entries of fii and f2 2 are real. If C is real-valued, that is to say all the coefficients in 
(|2.f p are real, pseudo-periodicity is equivalent to J\ = and —P~ X Q = fi 2 for some 
f2 G R dxd . In that case, the function z(t) may be returned as 
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z(t) = cos (ifi)zi + sin (tfl)z 2 — Q 1 J3, 
so that the second assumption in (|4.1|) reads as 

, , /cos ail sinafA . _ ,. ., 

Remark 4.1. The extension to the case P~ 1 Q = +Q 2 and Ji = is straightfor- 
ward and in this case the formula involves cosh(£f2) and sinh(£f2) in z(t). 

Remark 4.2. Periodicity of z(t) is obviously equivalent to exp(iTfifc) = 7^, for 
some T > and for all A; = 1, 2. 

Remark 4.3. The problem of existence of square or higher roots to real or complex 
matrices, as in (|4.ip . has led to huge bibliography. For instance, a simple criterion 
depending on elementary divisors for a real nonsingular matrix M to have real square 
roots is that each elementary divisor corresponding to a negative eigenvalue occurs 
an even number of times, see |5j pp. 413, Theorem 5]. But this result has been 
improved by Higham, since he proved that at most 2 r+c real square roots of a real 
nonsingular matrix M may be expressed as some polynomial in M ([51 pp. 416, 
Theorem 7]), r (resp. c) being the number of real (resp. distinct complex conjugate 
pair of) eigenvalues of M. 

4.2. Generation of pseudo-periodic solutions of D.E.L.. Let us study the 
existence, unicity and pseudo-periodicity of the solutions of (|2.4j) . In order to express 
the components of the solution of D.E.L. as in (|4.2|) . we use the main results in Section 
3 by adding the assumption that C is stationary. In that case, for all n € Is defined 
in (|3.3p . the matrix A n and the vector b„ do not depend on n. We set A = A n and 
b = b„ for n £ Is- 

Proposition 4.2. We suppose that 

\Sp(A)\ = ANd, l$Sp(A), M^^—^eN, 
and A3.10\) holds. Then the restriction of any solution x e to Q a ^ £ PI [a + 2Ne, b — 2Ne] 
is uniquely determined and its components have the shape 

g £ (t) = cst+ ^2 X! cstfcj,e(d fc )j exp ( — (t-a)j . (4.5) 

k = a,b c X p(i0)eSp(A) ^ £ ' 

je[l,cfl 

Moreover, the restriction of x e on Q a e l~l [a + 2Ne, b — 2Ne] is pseudo-periodic if and 
only if Sp(A) C U. 

Proof. Let us define the two vectors J5 and Jq in C d by : 

4N 

J 5 = p- 1 (^D-e(l) + J 3 ), Je = (Id-Y,B i , n )- 1 J 5 . (4.6) 

cnc-n r-. \ 

i—i 
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Note that J 6 is well-defined since 1 ^ Sp(A). When n G Ig, formula (|2.4[) may be 
rewritten under the form 

AN 

x(i + ne) = Si,„x(t + (n - i)e) + Jb- (4-7) 

i=l 

A particular constant solution of (|4.7[> is obviously given by x(i + he) = J6- We 
now apply Proposition 13.21 Given w € C d and v g c 4Jv ^ the vector sequence (A'w); 
satisfies the homogeneous recurrence (14. 7[) if and only if A € Sp(^4) is associated to v. 
Since A is diagonalizable, the eigenvectors are linearly independent and we get 

x £ (t +ne) = J2 X 7 w i + J 6, (4.8) 

where Wj are appropriate vectors in C d . If £ = (6 — a)/M, d s is a well-defined linear 
combination of d a and d;,, as seen in (|3.9p . Let us introduce the linear system of ANd 
equations 

AN d 

} j A™Wj = x e (a + Tie) - J 6 
j=i 

where the r.h.s. are computed from (|3.8|) . The determinant of this system is the Van- 
dermonde V{\\, . . . , A4at,2) which is nonzero since the eigenvalues {Xi}i of A are pair- 
wise distinct. Hence, due to p.9[) , the vectors wi, . . . , w^Nd are well-defined and may 
be uniquely written as linear combinations of d Q and dj,. If we denote the eigenvalues 
of A by A = exp(i8) with 9 € C, (|4.8j) may be rewritten as (|4.5|) . Pseudo-periodicity is 
equivalent to the requirement that 9 £ K for all exp(i#) g that is Sp(A) CU. I 

Proposition 4.3. Under the assumptions of Proposition \4~J^ and the hypothesis 
Sp(A) C U and e < (b— <z)/(4iV(d+l)), we may associate to any solutions : [a, b] — > 
C d of D.E.L. one and only one function y £ : [a, b] — >■ C d such that 

• y £ is a solution of D.E.L. on [a, b], 

• y e is pseudo-periodic on [a,b], 

• x e and y £ agree on Q a . £ H [a + 2iV£, 6 — 2iVe]. 

// x £ is pseudo-periodic then y £ = x,. : . Moreover, if x £ is continuous on [a,b], then, 
for all S > 0, sup tg [ a+(5 b _,5] ||x E (t) — y e (i)|| tends to as e tends to 0. 

Proof. Indeed, y e is generated by using (14.51) outside the grid and outside [a + 
2Ne, b - 2Ne], so that obviously x e = y £ on Q a>s n [a + 2Ne, b - 2Ne}. It turns out 
that y £ is also a solution of D.E.L. since the coefficients of the recurrence in (|2.4p 
are independent on time, that is to say the coefficients are the same for any grid. 
Due to the assumption Sp(A) C U, y £ is pseudo-periodic. Let us prove the unicity : 
we assume that there exists two pseudo-periodic solutions y £j i and y E ,2- Let us fix 
k € [l,d]. The component of index k of y e .2 — Ye l is of the shape (|4.5I) . So we 
may define 8 P € C as the coefficient of exp(i9 p (t — a)/e) in y E ,2(i) — ye,i(£) f° r au 
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p e [l,4Nd\. Suppose now that e < (b - a)/(4N(d + 1)). Setting t = a + (2N + n)e 
in (|4.5I) with 1 < n < ANd we get a linear system of size ANd such as 

4Af<i 

4 exp(m6» fe ) = (y e , 2 (a + (27V + n)e) - y M (<z + (2N + n)e))j = 0. 

fe=l 

By assumption, the Vandermonde determinant of this system is nonzero and we get 
S p = for all p. Since this holds for all component of y e ,2(t) — y £ ,i(t), we get unicity 
that is y e ,i(t) — y e ,2(t) for all t £ [a,b]. As a consequence of unicity, if x e is itself 
pseudo-periodic, then y £ = x e . 

Finally, let us choose e so that 2Ns < 8. Since x E and y E are uniformly continu- 
ous on [a + 8, b — 8], we choose s less than a modulus of uniform continuity for 8/2. 
If t E [a + 8, b — 5} and to = a + ne is the closest point of the grid to t, the triangle 
inequality yields ||x e (f) - y B (t)|| < ||x £ (t) - x £ (i G )|| + ||y E (te) - y e (t)|| < 8, ■ 

Remark 4.4. If the coefficients Ci are chosen as ji/e then the matrix A is a 

quadratic polynomial w.r.t. s. The eigenvalues of A are algebraic functions of e. 

Determining if Sp(A) is included in U is a polynomial elimination problem. For 

instance, if N = d = 1, the operators D e for which the spectrum of A is included in 

[— -1 [i —ll 

the unit circle are of the shape D e = O e 2 ' 2 + ikD E ' where k € M, see [6j pp.7, 
Proposition 5.2]. 

5. Obstructions to convergence of y £ (t) to z(t) as e tends to 0. 

5.1. Preliminary discussion. Under the assumptions of the three propositions 
of the previous section, to prove that 

y £ (t) tends to z(t) uniformly locally on ]a,b[ as e tends to 0, 
is not an easy task. It relies on the comparison of the formulas (|4.2[) and (14. 5|) . This 
is why we focus on phases and amplitudes occuring in z(t) and y £ (t). 

The convergence of y £ (t) to z(t) as e tends to is related to the three following 
properties. 

(a) If Xj is an eigenvalue of A which tends to 1 as e tends to 0, its phase 6j is 
such that -i- tends to a phase lo^ of some eigenvalue exp(iwfe) of SI. 

(b) For any phase weK such that exp(zo;) € Sp(£l), let f u £ C be the amplitude 
of exp(iuj(t — a)) in (|4.2|) . Similarly, let g £<U j £ C be the sum of the amplitudes 
occuring in (|4.5[) corresponding to exp(i6>) £ Sp(A) with 0/e ~ lj as e — > 0. 
Then lim |p e u - f u \ = 0. 

e— >0 

(c) The sum of the contribution gg in (|4.5|) of eigenvalues exp(z6') £ Sp(A) not 
tending to 1 cancels, as e tends to 0. 

Summing all triangle inequalities 



ge,uj exp ( i9- — - ) - f u exp (iu(t - a)) 
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UJ 

£ 



(b-a) 



over the group of eigenvalues tending to 1, and considering the contribution of eigen- 
values which are not tending to 1, ||y e (t) — z^Wc^ ([ a +s,b-S]) is upper bounded by 



£ 

' £ Sp(A) 

e /> o 



+ £ \9e*,-f u \+Me J2 I/- 

e™es P (n) 



6 Sp(A) 

e Sp(fi) 



sup 

~ eu -y 



U) 



If the three properties hold, the previous bound tends to as e tends to 0. Note lastly, 
that the result of convergence itself is related to the success of the shooting method 
and the convergence of the scheme. 

We shall illustrate in the following two subsections the convergence issue by giving 
two convenient examples when N — 1 and d = 2 for two special cases of □i 1 *' 15 '. In 
that case, we denote x„ = x(a + ns) = (x n , y n ). 



5.2. First example. Let us consider r £ 



restrict ourselves to the case where P = i^ 1 


P2^ 


■«-(s 


92 \ 




Pi) 




Qi] 


condition (|4.1[) implies that the two numbers 


qi - 


Vqi , qi 
and — 


- 92 


Pi - 


VP2 Pi 


~P2 



and Df ' r ' defined as in (|2.5[) . We 
and J x = 0. The 

are negative. In 



that case, we find that Sp(fl) = 



qi + q2 



Pi +P2 



qi - 92 



Pi -P2 



The recurrence v n +i 
note that 



Av n + b splits into two recurrences for X2 n and X2 



n+l- 



We 




j det(exp(if7(6 — a)) — exp(— i£l(b — a)))\ 

' ' ' ' ' ' ' P2 

(5.1) 

and accordingly to Proposition ^. 11 the shooting method is successful for z(t) if and 
only if the two eigenvalues of fi are not commensurable with tt. We get so far 



X n +2 



Vn+2 



= 2 



£ 2 Piqi - P2g2 \ _ , e 2 Piq2 - P2qi 



2 
Pi 



Pi 



Pi 



p\ 



U n •E n — 2 i 



£ Pi 92 - P29i 



The coefficients occuring in the previous recurrence are the entries of block i?2,n 
defined as in (|3.5p . Note by the way that the two blocks Bi n and E>3 n are zero. 
The sequences ((x2n> V2n)) and ((x2n+i,J/2n+i)) obey to the same recurrence but are 
computed independently each to the other. If M is even, the Dirichlet conditions for 
n = and n = M ensure existence and unicity of ((x2n,y2n)) provided the shooting 
method is successful. By reordering the components of the vector v„, the matrix A 

Q4 K J where K A = ^ Q 
Ok is the zero matrix of size k. Now the spectrum of K4 consists of the four numbers 

exp(±z0i), exp(±i6 l 2) where 

.2 



2 
Pi 



e 1 Piqi 



Pi 



P2q2\ 

— I Vn ~ Vn-2- 



pI 



is equivalent to a block diagonal matrix 



and 



?i = arccos 1 

2r 2 



qi + 92 



Pi +P2 



and 62 — arccos 1 

\ 2r 2 
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9i - 92 



Pi -P2 



We have here = 4 and 4Nd = 8, so Proposition (I4.2[) does not apply and 

indeed, the sequence ((x2n+i, 2/2n+i)) is not uniquely determined. Lastly, we get 

e 




qi ± qi 



Pi ±P2 



27 W1 < 2 - 



We note that property (a) holds if and only if r = — . The property (b) is much more 
delicate and is discussed in the last section. At last, property (c) is obviously true. 

5.3. Second example. We consider the operator used by Cresson in [3] to define 
scale derivatives : 

□e 2 ' 2 1 = -xi(t)^— x(t - e) + -x i + x-i i — — x * + £)• 

The characteristic polynomial of A may be factored into two biquadratic equations. 
The eight eigenvalues of A may be written as 



+ <2irV 1 - (ew fc ) 2 - CiV 1 + 2(ew fc ) 2 , 



2 ' "~ 2 

where (k — 1,2) are the eigenvalues of the matrix fi and Ci = Cf = 1- We see that 
the eigenvalues of A G C 8x8 are all distinct and of modulus 1. Looking for the limits 
as e tends to of the eigenvalues, we get four limits equal to 1, two equal to i and 
two equal to — i. The first four eigenvalues check the property (a) as shows expansion 
with Taylor series w.r.t. ewfc. Note that the four eigenvalues tending to 1 (obtained 
by choosing d = 1) may be written as A = exp(j^2Wfc) + 

The sub-sum of the terms in (|4.5|) implying eigenvalues which tend to 1 as e tends to 
may be rewritten as 

cskA 1 + l ^ £ C - \u 2 jZ 2 + 0(e 3 ))^, 

j=l,2,C=±l 

where 0(e 3 ) is uniform in t, by using Puiseux expansion of each factor of det(A— 
around e — 0. The limit of this sum as e tends to 0, is a combination of exp(iw / t(t — a)) 
and exp(— — a)) which is a step towards property (b). 

[id i+ii 

In fact, if we want to justify the choice of D £ 2 ' 2 , we may generalize a little bit 
the previous calculations to D e such that the operator of the l.h.s. of (|2.4p converges 
to (12. 3|) and is such that Sp(A) C U. In 6, Proposition 5.2], we prove that D e must 
be of the shape D^' s ^ with r = i — it and s = ^ + it. A formal computation of the 

eigenvalues Xj € Sp(A) shows that for two indices j\ and two indices j% we have 

it it 
lim A,, (e) = 1 + sqn(l — 2t)-, and lim A,,(e) = 1 — sqn(\ — 2t)-. -. 

e^O ■ /IW " y \+t \ + t 

[i ii tl=i i+ii 

The two operators D £ 2 ' 2 and D £ 2 ' 2 are the only ones such that 

N = 1, □e(l) = 0, n e (t) = 1 in [a + 2Ne, b - 2e] 

and Sp(A) C U for all C. As a conclusion, property (a) of Subsection 15. II holds. 
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6. Convergence of solutions and no n- resonance. This last section is de- 
voted to the convergence of y s (t) to z(t) in the case of multidimensional harmonic 
oscillator with N = 1. For sake of conciseness, we suppose that 

• each lagrangian is real, stationary, with Jk — 0, for all fc, that is to say 
£(x, x) = 5-Px 2 + |Qx 2 , with P and Q real, constant and nonsingular, 

• the operator D e is of the shape \3e with rel. 

The setting is discussed at the end of the section. Even with those restrictions, the 
convergence is not unconditional w.r.t. r. 

Lemma 6.1. Let us assume that limy £ (i) = z(i) for all lagrangian such that the 

e -»0 

hypotheses of the Propositions \4-l] and\4-2\ hold. Then, we have for some r, s G C 



□ £ = nl r ' 5] . 

Proof. Indeed, y £ (t) and z(t) exist, are unique and pseudo-periodic on [a,b]. Let 
Co be the lagrangian deduced from L by removing the terms J2 and J3 by 0. Since 
the solutions y e and y £j o of D.E.L. associated to C and £0 tend to the solutions z and 
Zo of the respective C.E.L. we have y e — y £ ,o = — Q _1 (D- e J2 + ^3) which tends to 
z — zo = — Q^ 1 J3. So, we obtain lim Q" 1 J2 x D- e l — for all J2. As a consequence, 

□ _ e l = Co + c\ + c_i = 0. Thus, D e = □| r ' s ' where r — c± and s — —c-\. ■ 

Lemma 6.2. Le£ </?(e, £) 6e a convergent Taylor series in some polydisc ofC 2 such 
that <p(Q, 0) = 0. If 5 > and B(Q, S) C C is the disc of radius S, let g : B(0, 5) -> C be 
a continuous function such that g{e) ~ J /or some r G C* and |.g(e) — j\ is bounded. 
Then, for all entire function ip(C) in C 7 we have 

limi>(g(sMe,en))=i> (r^(0,0)I d +T^(0,0)o\ . (6.1) 

uniformly in any compact subset ofC dxd . 

Proof. Let |.| be a norm of algebra over the Banach algebra C dxd . We denote 
by <p(e, C) = Yli j a i,j e% ^ tne Taylor series of <p(e, () in the bidisc |e| + |C| < r of C 2 . 
Let B'(0, n) G C dxd be an open ball of radius n > 0. We set Si = | min f 5, j^-) ■ 
The matrix-valued mapping (e,f2) i-> 8 := <y9(e,ef2) is well-defined and analytic in 
B(0,6i) x S'(0,ri). We have by assumption a ,o = and we denote ai,o = §f (0,0) 
and a ,i = §f (0,0). We get 

||ff(e)e-(TOi,o/,i+rao 1 in)|| < | ff (e) - j|||6|| + |r| ^ {aijWe]^- 1 ]]^. (6.2) 

-i+j>2 

Let 61 > such that \g(e) — -|| < 61 in [— S, 5]. Let JC be some compact subset of 
B(0, Si) x B'(0, ri). We use the fact that, in any bidisc |e| + |C| < r 2 < r , the Taylor 
series J2i j \ a ij\ £t C 3 i s normally convergent. The function Q/e — ip(e,eCl)/£ being a 
Taylor series, we have ||0|| < b2\e\ in K, for some 62 > 0. Similarly, the sum in (|6.2j) 
is upper bounded by £>3 1 £• | in K, for some 63 > 0. Then, the l.h.s. of (|6.2[) is upper 
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bounded by e(6i£>2 + l r l&3) in This implies that the convergence in (|6.1I) holds and 
is uniform in the ball B(0,ri) provided we have ip(C) = C- Now, if "0(C) is an entire 
function, the l.h.s. (|6.1[) gets sense and it is classical analysis that composition limit 
law holds for uniform convergence. ■ 

Theorem 6.3. We assume that M = (b — a)/e G N* is even. Let us consider 
□ e = Dg , for some r G R*. Let C be a real stationary quadratic lagrangian of the 
shape £(x, x) = \P~x? + ^Qx 2 , with P and Q real, constant and nonsingular. We 
suppose that —P~ 1 Q = f2 2 and Sp(fl) C U for some matrix f2 G R dxd diagonalizable 
over R. We require also the non-resonance property 

1 1 / (b — a) 2 ui 2 \ 

-(b-a)u(£Q and - arccos 1 - )4 Q, Vn € N*. (6.3) 

7r 7T \ 2r z n z J 

for all eigenvalue exp(iw) m Sp(f2)nU. TTien y E (i) is well-defined. Moreover, r = ±^ 
if and on/?/ if for all d a , db G y £ (i) tends to z(t) uniformly on [a, 6] as e — > 0. 

Proof. The necessary conditions of the first order C.E.L. and D.E.L., that is (|2 .3[) 
and (|2.4p . simplify into 

-Pz(t) + Qz(t) = 0, PO_ e n e y £ (t) + QYs(t) = 0, (6.4) 

completed with Dirichlet conditions. Let z(t) be the solution of (|6.4j) . that is to say 

z(i) = cxp(if2(t — a))fi + cxp(— ifl(t — a))f 2 , 
where t G [a, 6] and il 2 = —P~ 1 Q. Due to (|6.3[) . the diagonalizable matrix sin((6 — 
a)Q) is invertible and we may solve the boundary conditions z(a) — d a and z(b) = db. 
If we set 

F(r,d a ,db) = ^(sin(rr!)- 1 (exp(- J rr!)d a - d b ) (6.5) 

we get fi = F(6 — a, d a , db) and f 2 = F(a — b, d a , db). 

Since we have c\ — r/e, c_i = — r/e and Co = 0, then, for all t G [a + 2e, b — 2e], 
D.E.L. in (|6.4[) may be simplified into 

x(t + 2e) + x(t - 2e) = 2 ^ + ^ P" 1 ^ x(t). (6.6) 

Since f2 is diagonalizable over R, for some £? G R dxd , we have B _1 f2B = diag(u>i). 
For £ > such that e < 2\r \ min |u?j|, the matrix 



e i e 2 



6 = Sdiag I arcsin I 1 - J e R dxd (6.7) 

is well-defined and the computation of cos(i? -1 8-B) gives cos(O) = Ld— ^Vt 2 . By 
setting u„ = x £ (a + 2ne), the equation (|6.6|) gives u„ + i = 2cos(0)u„ — u n _i. The 
solution of this recurrence is given by 
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x e (a + 2ns) — exp(in0)g! + exp(— mO)g 2 
for all n from 1 to hM — 1. 

To compute the vectors gi and g 2 we use the Dirichlet conditions for y £ and the 
invertibility of the matrices sin(fcO), k > 1. Indeed, it is ensured since we have 

det(sin(/c6)) = f[ sin ^ arccos - j 7^ for all k, M e N* 

help to (|6.3[) . Now, the equation (|3.2p gives for t — a and t = b the boundary 
conditions 

x £ (a + 2e) = d' = (j d - d a and x £ (& - 2e) = d' b = (/ d - ^O 2 ^ d b . 

Note that dj, and d b tend to d a and d b respectively as e tends to 0. Next, (13.2[) gives 
for i = a + e and t = b — e the vectors 

gl - I (sin ((±M - 2) 6))" 1 (exp (-* (§M - l) 8) d' a - ex P (- l 6)d b ) 

g 2 = i (sin ((|M - 2) 9)) _1 (exp(i9)d b - exp (i (ifcf - l) 6) dQ . 
We also need in the following 

gr = | (sin (iMO))" 1 (exp (-i±M6) d a - d 6 ) 

g 2 = i (sin (|M6)) _1 (d b - exp (i±MO) d a ) . 
Since the sequence u„ is well-determined, Proposition 14.31 and its proof show how 
x £ (i) may be extended to an unique continuous pseudo-periodic function y £ (t) over 
[a,b]. 

We will apply several times Lemma [6.21 with <p(e,C) = arcsin (~C\J 1 ~ IP'C 2 ) 1 see 
(|6.7p . Keeping notation of the proof of Lemma [6~2l we have a^o = and fflo,i = -< If 



g{e) = we get if}(g(e)@) -> ip(^Q,). Applying this result to ^(C) = cos(C) or 
sin(C) or exp(±i£) we readily obtain lim gj = lim gj = fj, j = 1, 2, where 

fi = F fc^da, d„) and f J = F d a , d b j . (6.8) 

Now, all these preliminaries being done, let us discuss the convergence of y e to z. 
We fix S > 0. If t £ [a + 5, b — 5], we choose <?(e) = [^r] to represent the integer n 
such that t = a + 2ne. Hence, with t = t — a, g(e) ~ ^ and |g(e) — ^| < 1. If we 
define Z by 

Z(n, r 2 , d Q , d 6 ) = (sin^fi))- 1 sin(rifi)db - (sin^ft))- 1 sin((r 1 - r 2 )0)d a (6.9) 

for suitable numbers ri and r 2 , then we have z{t) = Z(t — a, b — a, d a , db) and we are 
going to prove that y e (i) tends uniformly locally on ]a, b[ to 

z'(t):=z(^,^,d a ,d b ). (6.10) 

Let us note that z'(t) does not stand for the derivative that we have denoted z(t). 
The vector z'(t) — y e (t) may be written as 

exp(ir£l)f[ +exp(— irfi)^ — (exp(mO)gi + exp(— in0)g 2 ). 
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Let us introduce the quantities 



qi 

02 

q:s 



-exp(m9)(gi - gi) 
exp(-m6)(g 2 - g 2 ) 



2r 

T 



= I cos ( — 91 ) — cos(n9) ) d Q 



2r 



q 4 = isin — -9l) (f[ - f£) - ism(n9)(gi - g 2 ). 



Straightforward computations yield z' (i) — y e (i) = qi 
the two equalities {[ + i' 2 = gi + g 2 = d a . Let 



that lim q 

e-s-0 



q4, as a consequence of 
v\\ be any norm on C d . Let us prove 
uniformly in [a, b], for i = 1, . . . , 4. As we have seen in the discussion 
before (|6.8|) . the vectors gi — gi tend to and exp(±m9) is bounded since 6 is real 
so the vectors qi and q 2 tend to uniformly on [a, b]. The case of q3 is obvious by 
using Lemma 15721 while q4 is more complicated. We obtain that ||q4|| is less than 

\da\\ + 



(sin (^?n)) 1 cos (^H) sin (£fi) - (sin (§M0)) 1 cos (±M9) sin(n9) 



(sin(^fi)) ^in(^O) - (sm(±M9)) 1 sin(n9) .||d 6 || 



that is 



cotg 



b-c 
2r 



■9. 



cotg I \m& 



|d«l 



cosec 



b-c 
2r 



cosec 



-M9 



■l|d 6 ||. 
(6.11) 



We have 11911 : \\B\\\\B 



(1 — 472 w ?)) where B diagonalizes 9. So, 
Lemma (|6.2|) with V'(C) = C0S (C) an d next ip(0 = sm (C) shows that the bound (|6.11|) 
tends to for all t € [a, b]. But this convergence is also uniform in [a, 6] due to formula 
(16. 2[) . to the previous bound of 9 and to the boundedness of g(e) — Hence, by 
using notation in (I6.10[) . we have proved so far that 

y £ (t) — > o z'(t) uniformly in [a, b\. 
Let us show that, if r ^ ±5, then we can choose t £ [a, b] and the vectors d a , dj, g C d 
in such a way that z'(i) 7^ z(i). Indeed, inspection of (|6.9|) shows that the coefficient 
of d a is the matrix 

f sin f^T fi ) ) sin ~ (sin (( 5 ~ a ) fi )) _1 sin ((* - ( 6 - 12 ) 

where sin (^r^) and sin((6 — a)fl) are invertible due to (|6.3[) . We choose i such that 
for all W € Sp(9) we have 

(sin (^u;))" 1 sin (^w) - (sin((6 - a)^))- 1 sin((i - 6)«) ^ 0. 
In this way, given any vector d^ we may choose d a so that z'(t) — z(t) = e-y where e\ 
is the first vector of the canonical basis of C d since the matrix (]6. 12[) is invertible. 
If r = ±1/2, we have z'(t) = z(t) and the proof is complete. ■ 



Remark 6.1. If M = (6 — a)/e is odd, the system (|3.9|) for determining d s = 
x e (a + e) from d a and d& is not Cramer. Indeed, D.E.L. in (|6.4[) simplifies into (|6.6[) . 
Since this recurrence does not match a, a + e and 6, the matrix occuring in (|3.10[) is 
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the zero matrix Od- So, when e is not of the shape alluded in the previous theorem, 
the convergence is not guaranteed (see the numerical experiments below). 

Remark 6.2. The assumption = 0, k > 2, is not restrictive. Let C be defined 
as in (|2.1|) and £o(x, x) = 5-Px 2 + ^Qx 2 . Let y e , y £! o, z and z be the solutions of 
D.E.L. and C.E.L. for £ and Co- We shall see that the property of convergence holds 
for C iff it holds for Co- Indeed, we have the formulas y e = y e ,o — Q~ 1 {0- E J2 + J3), 
z = zo — Q _1 J3 and lim y £ = z o- So we have lim y £ = z lim y e = z o- 

e->0 ' e->0 e->0 

Remark 6.3. We already obtained in [6] two characterizations of D^'^ among all 
operators D e . The first one was linked to the convergence of the l.h.s. of (|2.4j) to the 
l.h.s. of (|2.3p for all lagrangian (|2.ip and led to the relation r + s = 1. The second 
one ensured that Sp(A) C U if d = 1, which leads to s = r. Now, we may add a third 
family which consists in operators ni r,s ', for which the five-terms recurrence (|2.4[) 
splits into two three-terms recurrences, one for X2 n and X2„+i, this being equivalent 
to r = s. 

7. Numerical experiments. Let us illustrate the phenomenon of convergence 
proved in Theorem 16. 3[ when M increases. We set in every example below d a = 12, 
d b = -14, p = 1 and q = -0.23 

Figures [7.11 below illustrate the behaviour of y e (t) and z(t) when M increases 
(M = 30 and M = 120). We choose first 71 = 1/2, a = and b = 30, so condition 
(|6.3|) is true. 




Fig. 7.1. Convergence in non-resonant case (M = 30 and M = 120 ^ 



As soon as the condition (|6.3|) for the continuous lagrangian fails, the convergence 
does not occur. However, if b = a + ^(arcsin(p) + 2Kir), for p tending to 0, the 
upper bound (I6.11|) grows to infinity. The small denominators sin(M0) and sin(rw) 
occuring in y e (t) and z(t) imply that the convergence holds but is slowed down (see 
Figures 

Next, let us present two examples of non-convergence of solutions, that is to say 
when r = 71 7^ ±i. We remind that the study of the convergence of the operators 
in D.E.L. to the operators in C.E.L. is studied in [SJ pp.7, Theorem 6.3], and may be 
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Fig. 7.2. Slow convergence in quasi-resonant case (M = 1000 and M = 50000,) 



shortened as : □ = . Next, the existence of pseudo-periodic solutions implies 
r + s = 1 (see (i pp.7, Proposition 5.2]). If r = s = 71 ^ h, then □[ r ' r ' does not 
fullfil the previous requirements. We display in Figures 17.31 two such instances with 
71 = 0.6 and 71 = (l + i)/2. 




Fig. 7.3. Non-convergence phenomenon (71 = 0.6 and 71 = (1 J r i)/2) 



Let us conclude this paper with the following problems. First, formal and nu- 
merical codes have been written to do experiments on D.E.L. and C.E.L. in higher 
dimension (d > 2), to deal with huge matrices A and to work with non-periodic so- 
lutions. However, mainly due to the characteristic functions Xi(t) occuring in D e , no 
general pattern has been found neither for the convergence nor for non-convergence. 
Second, it would be interesting to get qualitative properties as continuity or mesura- 
bility of solutions x e (t) of D.E.L. as it is usual in the theory of functional equations. 
Another work is to relate the convergence of schemes to the convergence of solu- 
tions, none of these properties implying the other. These directions seem to be some 
interesting perspectives for subsequent work. 
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